Theory for the phase behaviour of a coUoidal fluid with competing interactions 
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We study the phase behaviour of a fluid composed of particles which interact via a pair potential 
that is repulsive for large inter-particle distances, is attractive at intermediate distances and is 
strongly repulsive at short distances (the particles have a hard core). As well as exhibiting gas- 
liquid phase separation, this system also exhibits phase transitions from the uniform fluid phases 
to modulated inhomogeneous fluid phases. Starting from a microscopic density functional theory, 
we develop an order parameter theory for the phase transition in order to examine in detail the 
phase behaviour. The amplitude of the density modulations is the order parameter in our theory. 
The theory predicts that the phase transition from the uniform to the modulated fluid phase can be 
either first order or second order (continuous). The phase diagram exhibits two tricritical points, 
joined to one another by the line of second order transitions. 
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I. INTRODUCTION 

The interactions between colloidal particles can be 
complex and varied. In order to understand and de- 
termine the interactions one must not only consider the 
direct or 'bare' interactions between the colloids them- 
selves, but also the influence of the medium in which they 
are dispersed (the solvent). For example, if the colloids 
carry a net charge g, they will tend to repel each other. 
However, the strength of this repulsive interaction is de- 
termined not only by the magnitude of the charge g, but 
also by the concentration and types of counterions in the 
solvent. The counterions may condense on the surface of 
the colloids, forming an oppositely charged double-layer 
around the colloids, which screens the Coulomb inter- 
action. Thus the range and strength of this repulsive 
interaction may be determined by controlling the concen- 
tration and type of colons and counterions in the solvent 

m. 

The interaction between the colloids may also be influ- 
enced by the presence of other species in the solvent, 
such as polymeric macromolecules. If these polymers 
adhere to the surface of the colloids they form a soft 
layer surrounding the particles thus effectively increas- 
ing their size. Such adsorbing polymers may be used to 
stabilise the colloids in suspension. On the other hand, 
non-adsorbing polymers have the effect of generating an 
effective attraction between the colloids. This depletion 
interaction arises due to the fact that when the sepa- 
ration between two colloids becomes less than ^ 2i?g, 
where Rg is the polymer radius of gyration, then the 
polymer chains are unlikely to be found between the two 
colloids since such confinement entails an cntropic cost. 



This depletion results in an unbalanced osmotic pressure 
on the colloids by the polymers, generating an effective 
attraction between the pair of colloids [J, [J. 

In the last few years there have been a number of stud- 
ies of systems of charged colloidal particles dispersed in 



a solvent containing nonadsorbing polymers 
Due to that fact that in these systems the charge is only 
weakly screened, there is a competition between a short 
ranged (depletion) attraction due to the presence of the 
polymers in the solution and a longer ranged (screened 
Coulomb) repulsive interaction. The competing interac- 
tions give rise to microphase separated fluid phases - i.e. 
to fluid states exhibiting equilibrium modulated struc- 
tures such as clusters or stripes. Other mechanisms lead- 
ing to this kind of competition are also possible. For 
instance, the attraction may be due to dispersion forces 
instead of osmotic depletion. Even with neutral particles 
interacting solely via depletion forces, competition may 
arise as a consequence of mutual interactions between 
the depletants, which modify the interaction potential 
with respect to the simplest Asakura-Oosawa picture 
Structuring and pattern formation has also been seen in 
two dimensional systems 0, where the repulsion is 
likely to be due to dipole-dipole interactions between the 
particles at the surface. Together with these experimen- 
tal studies, there has been signiflcant theoretical interest 
in model colloidal systems exhibiting such competing in- 
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An overview of the 
results and conclusions drawn from this body of work is 
given in Ref. [l^], and will not be repeated here. 

A simple model for these colloidal systems was consid- 
ered by some of us previously d, [13, El [13] ■ We model 
the 'bare' interaction between the colloids as a hard- 
sphere interaction, the long ranged screened Coulomb 
repulsion with a Yukawa potential and we model the de- 
pletion attraction using an additional Yukawa potential. 
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Thus the pair potential in our model is as follows: 

v{r) = Vhsir) + w{r) (1) 
where Vhsif) is the hard-sphere pair potential: 



Vhs{r) = 







r < a 
r > a, 



where a is the hard sphere-diameter and 

w{r) = —ecr h Act 



(2) 



(3) 



The first term in w(r) is the (attractive) depletion in- 
teraction and the second term is the (repulsive) screened 
Coulomb interaction. The magnitudes of the parame- 
ters e and A are determined by the concentration of non- 
adsorbing polymers in the solution and the charge carried 
by the colloids, respectively. The parameter Ai ~ cr/Rg 
is roughly the size ratio between the colloids and the 
polymers and the parameter A2 is determined by the con- 
centration and type of screening ions in the solvent. In 
the present work we consider the case when Ai > A2, so 
that the pair potential v{r) is repulsive for large r [32l |. 
However, at intermediate distances r > a between the 
particles, v{r) may be attractive - i.e. v{a^) < 0. This 
occurs when the parameter e > A. In this case there is 
a competition between the long range repulsion and the 
shorter ranged attraction. Depending on the values of 
the set of pair potential parameters {e. A, Ai, A2}, the 
fluid temperature T and number density p, this model 
may exhibit various modulated phases [ll|, [ij] . 

The precise nature of the phase transition from the uni- 
form fluid phase to the modulated fluid phases is not well 
understood. In Ref. [ll| some of us developed a density 
functional theory (DFT) [H, |33| for the present model. 
It was found that the DFT provided a good account of 
the structure and thermodynamics of the uniform liquid. 
The theory also predicts that for certain choices of the 
parameters {e, A, Ai, A2} and the temperature T, there 
exists an interval of the fluid density p where the uniform 
fluid becomes unstable against periodic density fluctua- 
tions of a certain wavelength, indicating that there must 
be a phase transition to an inhomogeneous (modulated) 
phase. The locus of this line in the phase diagram is 
called the A-line 0, [H, [H, [33,j38[ . 

Following the study in Ref. [ll|, a further study was 
made using Monte-Carlo (MC) computer simulations and 
integral equation theories T^, in order to determine the 
nature of phase transition occurring at the A-line. The 
results from the integral equation theories were mostly 
inconclusive - all the closure relations that were imple- 
mented either failed to describe the phase transition or 
had no solution in the region of the phase diagram where 
the inhomogeneous phases are expected, with the excep- 
tion of the Percus-Yevick (PY) closure, which is able to 
describe to some extent the nature of the transition from 
the low density homogeneous phase to the modulated 



(cluster) phase. However, this theory also failed to have 
a solution at higher densities. 

The picture that emerged from the MC simulations 
was that for certain values of A and there is a first 
order phase transition from a low density homogeneous 
phase to a modulated (cluster) phase and also at higher 
densities a first order phase transition from a modulated 
(bubble) phase to the uniform liquid (l^ . It was found 
that there is a difference in the densities Ap between 
these coexisting phases. However, Ap was found to be 
a fairly small quantity, which means that the density p 
is not a good order parameter for this phase transition. 
We show below that the amplitude of the density modu- 
lations A proves to be a better order parameter for this 
phase transition. On following either of these two first 
order transitions to larger values of e~^, one finds that 
Ap ^ 0. The two critical points, at which these two 
transitions cease to be first order (i.e. when Ap = 0) are 
difficult to locate using MC simulations ■ One reason 
for this was due to a lack of knowledge of the relevant or- 
der parameter for the phase transition. Based on a com- 
parison with the phase behaviour of lattice models with 
competing interactions [H, [s^, it was concluded that 
the two critical points were most likely to be tricritical 
points [l^l , connected to one-another by a line of second 
order (continuous) transitions. However, it was not pos- 
sible to see any signature of the second order transition 
line in the MC simulation results, leaving the conclusions 
somewhat tentative. 

In the present work we use DFT to examine the phase 
behaviour in the vicinity of the A-line and we find that 
the theory confirms the scenario proposed in Ref. [l2] - 
i.e. that the critical points are indeed tricritical points 
and that connecting these is the A-line itself - a line of 
second order transitions. We also use the DFT theory 
to develop an order parameter (Landau) theory for the 
phase transition, by making a sinusoidal approximation 
for the density modulations and then expanding the free 
energy in powers of the order parameter A, allowing us 
to investigate the transition in detail. We also confirm 
our results numerically, solving the full DFT theory for 
the density profiles. 

This paper proceeds as follows: In Sec. |TT]we describe 
the simple mean field DFT that we use to study the phase 
transitions. In Sec. IIIII we determine structural proper- 
ties predicted for the uniform fluid phase, focusing in 
particular on the static structure factor. We find a sim- 
ple expression for the wave number at which the static 
structure factor has a peak, corresponding to the typical 
length scale of the fluid modulations. In Sec. lIVi starting 
from our DFT we develop an order parameter theory for 
the phase transitions in this system, by expanding the 
free energy in powers of A. From our approximate free 
energy we determine the phase diagram of the system for 
particular choices of the pair potential parameters. In 
Sec.|V]we determine the phase behaviour by solving the 
full DFT, showing that our simplifled order parameter 
theory captures most of the important physics. Finally, 
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FIG. 1: Phase diagram in the versus density p plane 
obtained from the LDA DFT ([7]) for the set of pair potential 
parameters Ai = 1, A2 = 0.5 and A = 0.5. Within the region 
enclosed by the A-line the uniform fluid is unstable against 
periodic density fluctuations with wavevector fc = fee 7^ 0. 
The liquid-vapour coexistence curve (the binodal) and the 
spinodal are obtained from the same free energy functional 
for the uniform fluid. 



in Sec. I VII we draw our conclusions. 



II. MEAN-FIELD DFT FOR THE MODEL 

For systems of particles interacting via potentials of 
the form in Eq. ([1]), a common approach is to split the 
Helmholtz free energy into a contribution from the hard- 
sphere interactions between the particles - the 'reference' 
part - and a separate part due to the slowly varying 
tail of the potential w{r) [s^l- Taking this approach, in 
Ref. the following approximation for the intrinsic 

Helmholtz free energy functional for the system was pro- 
posed: 

+ \J^^J dr'p(r)p(r')w(r - r'). (4) 

where ^/fs°^[p(r)] is the Rosenfeld approximation [4ll.l4^. 
I43I I for the hard-sphere contribution to free energy. The 
Rosenfeld functional is a non-local weighted density func- 
tional, which also includes the exact ideal-gas contribu- 
tion to the free energy. The remaining contribution is 
a simple mean-field approximation for the contribution 
to the free energy from the longer ranged interactions 
between the particles. 

Note that in the mean field contribution to the free 
energy ^ [or one requires the value of the poten- 
tial w{r) for all values of r. However, given the hard 
sphere contribution to the pair potential, Eqs. JI]) - ([2]), 
the value of w{r) for < r < a should be irrelevant. 



Thus, the value of w{r) for < r < a used in the mean- 
field contribution to Eq. is effectively a free param- 
eter in the theory. In Ref. [11] it was found that the 
choice w{r) — w{a~^) — —e + A, for < r < a gave 
better agreement with SCOZA for the fluid phase, than 
extending down to < r < u the double Yukawa form 
of w(r) that is used for r > a. In the present study we 
use the truncated potential w{r) = w{a~^) = —e + A, for 
< r < a for all the DFT calculations. However, in Sees. 
mil and IIVI we do not truncate u;(r) in order to simplify 
the analytic calculations. 

The equilibrium one-body density profile p{r) is deter- 
mined by minimising the Grand potential functional 



(5) 



where p is the chemical potential and Vea;t(r) is the ex- 
ternal potential. The equilibrium density profile is the 
solution to the Euler-Lagrange equation 



m[p\_ 

6p{r) 



= 0. 



(6) 



If one requires the theory to be able to describe the oscil- 
latory (with wavelength ~ a) density profiles that occur 
when the fluid is subject to an external potential that 
varies strongly over short distances (such as, for exam- 
ple, the external potential due to the wall of the fluid 
container), then one must use a weighted-density DFT 
such as the Rosenfeld theory. However, if one is more 
interested in the large scale structures (stripes, clusters 
etc) that arise in the present system, due to the com- 
peting interactions in w(r), then one may simplify the 
above theory, by making a local density approximation 
(LDA) in the reference hard sphere functional. In this 
case one may assume that the intrinsic Helmholtz free 
energy functional of the system is given by the following 
mean-field approximation: 



T[pir)] 



drp(r)/(p(r)) 



-/dr 



dr'p{r)p{r')w{r ~ r') (7) 



where f{p) is the Helmholtz free energy per particle of 
a uniform fluid of hard-spheres with bulk density p. We 
use the Carnahan-Starling approximation (3^ : 



f3f{p) = Hrj) 



^(4 - Stj) 



(8) 



where 77 = irpa^ /6 is the packing fraction and /3 — l/ksT 
is the inverse temperature, which we will henceforth set 
to he (3 — 1. From Eqs. (O - ([5]) we obtain the follow- 
ing Euler-Lagrange equation for the equilibrium density 
profile: 



/(p(r)) + p(r)/'(p(r))+ / dr'p(r')i«(r-r') + 14.t(r) = 

(9) 
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FIG. 2: In the upper figure we display the density profile for 
the fiuid confined between two parallel hard walls, separated 
a distance L — 120a, in the case when the pair potential 
parameters are — 0.65, A = 0.5, Ai — 1, X2 = 0.5 and the 
average fiuid density is pa'^ = 0.25. In the lower figure, we 
display a magnification of the left hand portion of the density 
profile. The solid (red) line is the density profile calculated 
using the LDA reference hard sphere functional, Eq. ([7]). The 
dashed (black) line is the result obtained using the Rosenfeld 
functional for the hard-sphere reference functional. 



The free energy of the uniform bulk liquid may be ob- 
tained by setting the density profile p{r) = p, i.e. a con- 
stant, in either of Eqs. or ([7]). The resulting free en- 
ergies differ slightly between the two functionals, since 
for a uniform fluid the hard-sphere contribution from 
the Rosenfeld functional ^ is equivalent to that from 
the scaled particle (PY compressibility) equation of state 
[S^ l , whereas the hard-sphere contribution to the free en- 
ergy from ([7]) is that from the Carnahan-Starling equa- 
tion of state. However, for the densities of interest here, 
the difference between the two free energies is small. 

In Fig. [1] we display the uniform vapour-liquid coexis- 
tence curve (binodal), the spinodal and the A-line, for the 
set of pair potential parameters: Ai = 1, A2 = 0.5 and 
A — 0.5. Here, as in Refs. [ill, [l3> ^^e phase diagram 
has been mapped by changing the attraction strength e 
at fixed repulsion strength A. This can be achieved in a 
solution of charged colloidal particles by changing the de- 
pletant concentration at fixed temperature and salt con- 
centration. This choice implies that the region of e values 
where competition is important and microphase forma- 
tion is expected has both a lower and an upper bound: 
at low e, the interaction will be mostly or entirely repul- 
sive, the A-line is not met, and no transition will take 
place, except for those involving the occurrence of solid 
phases. At high e, the interaction will be mostly attrac- 
tive, and bulk liquid-vapour phase separation will take 
over. This is shown in Fig. [U where the liquid- vapour 
binodal moves out of the region bounded by the A-line 
for values of > 0.55. 

In order to demonstrate the reliability of the LDA DFT 



versus the more sophisticated Rosenfeld DFT, we calcu- 
late the fluid density proflles for the fluid confined be- 
tween two parallel planar hard walls separated a distance 
L = 120(7, where Vextir) = Vext{z) = for < z < L 
and Vext{z) — 00 otherwise. The external potential varies 
only in the Cartesian z-direction, and we assume that the 
fluid density profile also only varies in the z-direction, so 
that the Euler-Lagrange equation ^ becomes 

p = f{p{z))+p{z)f{p{z))+ f dz'p{z')<j,{z-z') + Vext{z), 



(10) 

where (l){z) — ^ Ax J dywlr). We solve this equation for 
the equilibrium density profile p{z) by discretising the 
density profile and using a simple iterative numerical pro- 
cedure. In Fig. [2] we display results for the case when the 
pair potential parameters are chosen to be = 0.65, 
A = 0.5, Ai = 1, A2 = 0.5 and the chemical potential 
is chosen so that the average fluid density in the slit is 
pa^ = 0.25. As one can see from Fig. [T]this corresponds 
to a state point that is well inside the A-line, where the 
fluid is strongly modulated. Comparing the two density 
profiles displayed in Fig.[2l we see that they are are fairly 
similar. The main difference occurs for values of z that 
are close to the confining walls. The Rosenfeld functional, 
which gives a better account of the hard sphere correla- 
tions, predicts that the fiuid profile contains oscillations 
with wavelength cr in the vicinity of the wall, due to 
packing of the spheres at the wall. However, the density 
profile obtained from LDA functional ([7]) does not have 
these small length scale modulations and only exhibits 
the larger length scale modulations that arise due to the 
competing interactions in w(r). 



III. STRUCTURE OF THE UNIFORM FLUID 

From previous studies of the present model fluid d, [l^, 
[HI, [13, E E lis ill il [11, [23, [11, il, we know that 
the static structure factor S{k) plays an important role 
in characterising the microphase structuring displayed by 
the system. One finds that there is a large peak in S{k) 
at /c = fcc <C 27r/cr, where Ic = 27r /kc is the length scale 
associated with the density modulations in the system. 
S{k) is given by the following expression 

where c(fc) is the Fourier transform of c(r), the bulk fiuid 
(Ornstein-Zernike) direct pair correlation function [s^ . 
This function may be obtained from the free energy func- 
tional via the following relation [33, 3^: 



c(r,r') = -f3 



5p{r)5p{r') 



where 



J'^d[p{v)] = ksT J drp(r)[ln(A3p(r)) - 1] 



(12) 



(13) 
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is the ideal-gas contribution to the free energy; A being 
the thermal wavelength. For the homogeneous bulk fluid 
where p{r) — p, we find c(r, r') = c(|r — r'|) = c{r). From 
Eqs. ([7]) and we obtain the following approximation 
for the pair direct correlation function 



c(r, r') [2/'(p(r)) + p(r)/"(p(r)) - ^ 

—(3w{r — r'), 



r') 

(14) 



where /' and /" are the first and second derivatives of 
/ with respect to p and 6{r) is the Dirac delta function. 
Thus the Fourier transform of this quantity evaluated for 
the bulk fluid is 



c(fc) = -13 



2f{p)+pf"{p) 



■f3w{k), (15) 



where w{k) is the Fourier transform of w{r) and the first 
term in the right hand side is equal to —i/{pXhf) + VPj 
Xhf being the isothermal compressibility of the hard- 
sphere fluid divided by that of the ideal gas. This treat- 
ment amounts to taking the expression 



c(fc) = Chs{k) - (3w{k) 



(16) 



given by the random phase approximation (RPA) and 
setting fc = in Chs{k), the Fourier transform of the 
hard-sphere direct pair correlation finction. In Ref. [ll| 
the uniform phase was studied by the full RPA, which 
is obtained from the free energy functional ^ adopted 
there in the homogeneous regime. If we define the fol- 
lowing parameters: ci = eexp(Ai) and €2 = Aexp(A2), 
then we may rewrite the pair potential w{r) as follows: 



exp[-Ai(r/g)] exp[-A2(r/cr)] 
w{r) = -ei ; h €2- 



r/a 



From this we obtain: 



w{k) 



r/a 



(Af + )fcV2) (Ai + fcV2)' 
Thus the static structure factor is given by 

S{k) ^ 



pm'{p)+pf"{p)+Hk)] 
1 



D{ky 



(17) 

(18) 

(19) 
(20) 



which defines the denominator function D[k). For the 
sets of parameters {Ai, A2, ei, £2} that we consider here, 
the structure factor given by Eq. (fT9|) exhibits a single 
peak, ai k = kc- This is the peak that characterises the 
cluster/stripe modulated structures in the system. This 
approximation for S{k) becomes unreliable for large wave 
numbers k ^ kc, as one should expect, given the delta- 
function approximation for the hard-sphere contribution 
to c(r), in Eq. (fU)) . 

One finds that in a certain portion of the phase dia- 
gram, the uniform fluid is unstable with respect to peri- 
odic density fluctuations - this occurs when S{k) —>■ 00 



at fc = fee- In other words, when D{k = kc) 0. The 
locus in the phase diagram at which D{k — kc) = is 
the A-line [11]. In Fig. [l] we display the A- line for the 
set of pair potential parameters: Ai = 1, A2 = 0.5 and 
A = 0.5. We see that this instability pre-empts the spin- 
odal instability which corresponds to S{k = 0) ^ 00; i.e. 
when D{k = 0) ^ 0. 

The wavelength of the periodic density modulations is 
Ic = 27r/fcc; this length scale is the typical distance from 
the peak of one density modulation to the next. Within 
the present theory, we are able to obtain a simple expres- 
sion for kc, as a function of the pair potential parameters. 
To do this we recall that kc is the value of k for which 
S{k) is maximum. This maximum in S{k) corresponds 
to a minimum in D{k), i.e. when 



dD{k) 
dk 



0. 



(21) 



Given that within our approximation for D{k), Eqs. ([T^ 
and ((20| . the only fc-dependence in D{k) enters in the 
function w{k), we find that this condition simplifies to 
the following: 



dw{k) 
dk 



SiTCia^k 

(A? + fc2a2)2 



8TTe2(T^k 

2^2 



fc2cr2) 



= 0. 



(22) 



There are two solutions to Eq. ([2211 . The first is fc = 0. 
This corresponds to a maximum in D(k) and a minimum 
in S{k). The second solution to Eq. (I22p is fc = fcc, where 



kc = 




(23) 



and where a = yei/ez. Note that this expression for kc 
does not depend of the fluid density p, indicating that the 
wavelength of the density modulations is independent of 
the fluid density. However, the amplitude of any density 
modulations does depend on the fluid density, since this 
quantity depends upon the height of the maximum in 
5(fc), which does depend on the density p. 

We are now in the position to obtain a relatively simple 
expression for the value of ei on the A-line, as a function 
of density and the other pair potential parameters. Re- 
call that the A-line corresponds to the line in the phase 
diagram where D{k = kc) — [ll|. From Eqs. (fTO|) and 
([20)1 we see that this condition is equivalent to the con- 
dition 

2f{p)+pr{p) + w{kc)^0. (24) 
Using our expression for kc in Eq. (I23|) . we obtain: 

4:Tra^e2{a ~ 1)^ 



w{kc) 



(25) 



Eqs. ([23]) and (|25p may be solved in a straightforward 
manner, by fixing the value of the density p and then 
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form 



FIG. 3: In the upper figure we display kc, the wave-number 
for the maximum in S{k), calculated along the A- line, for the 
case when A = 0.5, Ai = 1 and A2 = 0.5. The solid line is 
the result in Eq. (|23[) and the dashed line is the result from 
the RPA. Below this we display the A-line. The solid line is 
the result in Eq. (|26p and the dashed line is the result from 
solving the RPA. 



solving Eq. ((24|) for ei, to obtain: 



£1 - £2 1 



\Xl-\l){2f'{p)+pf"{p)) 



(26) 



This gives the value of ei on the A-line, as a function of 
the fluid density, p. 

In the lower panel of Fig. [3] we compare the result for 
the A-line in Eq. with the result obtained from the 
RPA, Eqs. (HU and IH]). RecaU that Eq. jMl) is obtained 
by taking the RPA expression for the pair direct corre- 
lation function, Eq. (|16p . and setting A: = in Chs{k). In 
the upper panel of Fig. [3] we display the value of kc calcu- 
lated along the A-line, obtained from Eq. ((23)l (solid line) 
and also the result from the RPA, Eqs. pTjl and 
(dashed line). We see fairly good agreement between the 
two sets of results, providing further confirmation that 
the the LDA free energy functional ^ is able to account 
for the correlations in the fluid in the portion of the phase 
diagram in the vicinity of the A-line. 



IV. ORDER PARAMETER THEORY 

In this section, starting from the DFT proposed in Sec. 
im we develop an order parameter theory for the phase 
transition from the uniform fluid to the modulated fluid 
phase. To develop the theory, we make the following 
assumptions: (i) we assume that the density profile varies 
only in one spatial direction, i.e. we assume p(r) — >• p{z). 
(ii) We assume the density modulations Sp{z) are of the 



p{z) = p + Sp{z) 

— p + Asm{kz), 



(27) 



where p is the average fluid density and the amplitude 
A is the order parameter for the transition. We should 
expect the wave number for the modulations in Eq. (j27p 
to be A; = fee- However, for the present we will not make 
any assumptions concerning the precise value of k in Eq. 
(|27p . other than to assume it is finite and non-zero. 

From our assumption that the density profile varies 
in only one Cartesian direction we may re-write the 
Helmholtz free energy ([7]) as follows: 



T[p{z)] 



J dzp{z)f{p{z)) 

L2 



dz / dz'p{z)p{z')(f>{z - z') (28) 



where the limits on the integrals go from ~L/2 to L/2 
(we implicitly take the thermodynamic limit L 00), 
= J dx J dy is the system cross sectional area and 



iz) 



dx J dy'w{r) 



Ai 



A, 



(29) 



where we have used Eq. \Yl\ to obtain the second line in 

Eq. jig). 

Substituting Eq. (^7)1 into (pS)) and making a Taylor 
expansion in the LDA hard-sphere part of the free energy, 
we obtain: 



T\p\ = J dzjp, 



1 dHpf) 




2 9p2 

Ip^wiO) + ^wik)5p' 



(30) 



w{k) 



dz Sp^ 



dzSp'^ + Oi / dz6p^) (31) 



where we have used the fact that all terms which involve 
an integral over an odd power of 5p are zero, due to 
the fact that Sp is a sinusoidal function of z. If we now 
assume that the system size L = mn/k, where m is an 
integer, then we obtain: 

^-^-^-^^A^-^^,A^-,OiA^) (32) 



where V = is the system volume, 

= pf{p) + ^(^> 



(33) 



is the free energy of the uniform fluid and B — 
(3[d\pf )/dp% = /3(4/"'(p)+p/""(p)). The higher order 



terms in the expansion (1321) are fairly straightforward to 
obtain, since they are merely terms arising from the Tay- 
lor expansion of the hard-sphere part of the free energy 
functional. For example, the next term of order 0{A^) 
in Eq. is CA^/6, where C = l3[d^{pf)/dp%/-iM. 

We may now use Eq. ((5^ to investigate the phase be- 
haviour in the vicinity of the A-line. For a given state 
point, the equilibrium value of the amplitude A is that 
which minimises the free energy ((5^ - i.e. the equilibrium 
value of the amplitude is the solution to the equation 



dA 



0, 



(34) 



which together with ([5^ gives us the following equation 
to be solved for A: 



R^A+^A^ + 0{A'') = Q 



(35) 



We first note that the coefficient _B > for all densities 
< p < G/ttct'^ and the same is true for all coefhcients of 
higher order terms in the expansion. Secondly, we note 
that for the uniform fluid, outside the A-linc, D{k) > for 
all values of k. Thus, outside the A-line, the minimum of 
the free energy ((32|) is when ^ = 0, as one should expect. 

On the A-line itself, D{k) > for all values of fc, with 
the equality being the case only for k = kc- However, 
inside the A-line, D{kc) < 0. Thus the coefficient of A^ 
in the free energy (|32p is negative for k = k^ meaning 
that inside the A-line the free energy is lower for some 
^ > 0, than when ^ = 0, indicating that the modulated 
fluid has a lower free energy than the uniform fluid. If 
we neglect the terms 0{A^) and higher in the expansion 
([32]). then from Eq. p5)) we find that the minimum of 
the free energy is when the amplitude 



A 



^\D{kc)\ 
pB 



1/2 



(36) 



Substituting this value of A into Eq. ([32|) we obtain 



I3T 



V 



(37) 



where we have neglected the contribution from terms of 
0{A^) and higher in Eq. (|37p . Since we have fairly simple 
expressions for all the terms in Eq. ([57)1 [see Eqs. (HH]), 
([20)l . ((25)) and ([33])], we now have a simple expression for 
the Helmholtz free energy of the modulated phase. 

Before using Eqs. (l33|) and (I37p for the free energy to 
examine the nature of the phase transition between the 
uniform and modulated fluid phases, we flrst recall that 
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FIG. 4: The Helmholtz free energy per particle a plotted as 
a function of v, the volume per particle, for the set of fluid 
pair potential parameters A = 0.5, Ai = 1 and A2 = 0.5. The 
solid line is the free energy for the uniform fluid, given by Eq. 
(|33p and the dashed line is the free energy for the modulated 
fluid, given by Eq. p7p . Fig. (a) is for the case when = 
0.83, above the tricritical point, and Fig. (b) is when = 
0.72, below the tricritical point. The dotted lines show the 
common tangent construction between the coexisting state 
points, which are marked with the symbol •. The points on 
the A-line, where the two free energies are equal are marked 
with the symbol x . 



for two phases to coexist, the temperature, pressure and 
chemical potential in the two phases must be equal. As 
we show in the Appendix, these conditions correspond 
geometrically to making a common tangent construction 
on the free energy per particle a{v) — J- /N plotted as 
a function of z; — V/N — the volume per particle 
{N is the total number of particles in the system). In 
Fig. m we display the free energy a as a function of v 
for two different values of e~^. In Fig. [3)Ja) we display 
the free energy for the set of pair potential parameters 



Ai = 1, A2 = 0.5, A = 0.5 and e'^ = 0.83. This curve is 
typical of the case when > e^^, where e^^ is the value 
of at the tricritical point. In this case we see that 
there is no common-tangent construction - i.e. the free 
energy is a convex function. Thus the phase transition 
between the uniform fluid and the modulated fluid occurs 
at the A-line and the transition is second order. In fact, 
on increasing the density p (decreasing w), we flnd there 
are two second order transitions. These two points are 
marked by the symbol • in Fig. [D^a), which are also two 
points on the A-line. The low density tricritical point is 
at e^^ = e^^ ~ 0.75 and the higher density one is at 
e^j,^ ~ 0.79. 

On decreasing e^^ below e^^^, one finds that the free 
energy is no longer convex and that the common tangent 
construction between the uniform fluid free energy and 
modulated fluid free energy can be made - see for ex- 
ample the results for e^^ = 0.72, displayed in Fig. [IJb). 
The common tangent construction lines are the two dot- 
ted lines. These lines join coexisting state points which 
are marked by the symbol • in Fig. [4jb) which lie at 
points either side of the A-line (marked x). The value 
of e"^ at which the free energy curve goes from being 
convex, to non-convex defines - this is the tricritical 
point floj . 

In Fig. Owe display the phase diagram resulting from 
performing the common tangent construction on the free 
energy for a range of values of e~^. We also display 
the phase diagram when the pair potential parameter 
A = 0.1. In calculating these phase diagrams we in- 
cluded the term of 0{A^) in the free energy, since trun- 
cating the expansion of the free energy in powers of the 
order parameter A is only really justified when A is small. 
Even including the term of 0{A^), one should not expect 
the theory to be reliable for determining the coexistence 
curve between the uniform fluid and the modulated fluid 
for < e^^, away from the tricritical point. Owing to 
this, the phase diagrams in Fig. [5] are only qualitatively 
correct for < e^^, where the truncated theory pre- 
dicts the free energy of the modulated phase to be lower 
than it really is. For the case A = 0.5 the theory predicts 
(incorrectly) that as is decreased, the two first order 
transition lines never intersect, as they do for A = 0.1. 
This is the reason why the theory (incorrectly) predicts 
that there is no triple point for the case when A — 0.5. 

A final question we wish to address in this section is: 
How does the amplitude A grow as a function of distance 
in the phase diagram from the A-line? Let us denote the 
fluid density on a the A-line itself, for a certain value of 
> e^^, as p\. Recall that D{kc, p\) = 0, which means 
that D{kc) ^ {p — p\) near the A-line. Combining this 
with Eq. (|36p . we flnd that inside and near to the A-line 
the amplitude 

A^ip^pxf'^ (38) 
and, of course, outside the A-line the amplitude ^ = 0. 




FIG. 5: Phase diagram in the versus density p plane, for 
the system with pair potential parameters Ai = 1, A2 = 0.5 
and (in the top diagram) A = 0.1 or (in the bottom diagram) 
A = 0.5, obtained using Eq. (|32p truncated after the term of 
OiA""). 



V. DFT RESULTS 

In this section we present results from numerically solv- 
ing the full DFT theory ((7]), by flnding the solution to 
the Euler-Lagrange equation (jlOp . without resorting to 
the sinusoidal approximation that was used to develop 
the order parameter theory in the previous section. 

To calculate the density profile that minimises the 
grand potential functional for a given point in the phase 
diagram, i.e. for given values of e^^ and chemical poten- 
tial /i, we solve the Euler-Lagrange equation pU|) using 
a simple iterative scheme. To determine the equilibrium 
density profile for the bulk system, i.e. when Vext{z) = 0, 
we solve for the profile using periodic boundary condi- 
tions. Note that in this case a uniform (constant) den- 
sity profile p{z) = p corresponds to a stationary curve 
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FIG. 6: Equilibrium density profile inside the modulated 
phase region, for = 0.65 and pa^ = 0.29. The same choice 
of the potential range and strength was made calculating the 
phase diagram in Fig. [T] where we see that this density pro- 
file corresponds to a point well inside the modulated phase 
region. 



value of the density profile to be p at a single point. 
This was done for the following reason: If a particular 
density profile p{z), such as that displayed in Fig. [HI is 
the equilibrium density profile, then so is p{z + Z), where 
I is any real number. Thus, any numerical optimisation 
algorithm will initially descend to a minimum of the free 
energy corresponding (say) to the density profile p{z). 
However, due to numerical errors the gradient of the free 
energy landscape will not exactly be zero for this profile 
and so the optimisation routine will then perform a very 
slow walk along the free energy valley, where points along 
the bottom of the valley correspond to different values of 
I. By fixing the density at a single point, one penalises 
such translations. 



A further issue to consider is that since we solve for 
p{z) on a finite grid of length L with periodic boundary 
conditions, this length L must be commensurate with 
the periodicity of the modulations in the density profile. 
Thus we must also minimise the free energy with respect 
to variations in L, or, equivalently, with respect to vari- 
ations in the spacing between density grid points. 
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FIG. 7: Phase diagram in the e" ~ p plane obtained from 
the RPA DFT for the parameters Ai = 1, A2 = 0.5 and A = 
0.5. The shaded regions are the two phase coexistence regions 
between the uniform and modulated fluid phases. They end 
at the tricritical points, where the transition ceases to be first 
order and becomes second order. Above the tricritical points, 
the transition line is the A-line. 



for the functional ([7]). In order to find the minimum of 
the free energy corresponding to an oscillatory density 
profile, such as that displayed in Fig. [6l we must break 
the symmetry by using a non-uniform initial guess for 
the density profile in our iterative solver. The final den- 
sity profile does not depend on the precise form of the 
initial guess; we used either a step function or the den- 
sity profile from a neighbouring state point as our initial 
guess in this study. We also found it useful to fix the 



Phase coexistence between the uniform fluid and the 
modulated fluid phases is determined by calculating the 
grand potential fl for fixed values of (which is equiva- 
lent to fixing the temperature) , whilst slowly varying the 
chemical potential p. Recall that two phases coexist if 
the chemical potential, pressure (recall SI = —PV, where 
P is the pressure) and temperatures are equal in the two 
phases. We display the resulting phase diagram in Fig. 
[T] The shaded region denotes the two phase coexistence 
region. This phase diagram should be compared with the 
lower phase diagram in Fig. O which was obtained using 
the order parameter theory of the previous section (re- 
call, however, that Fig. [5] was obtained using a different 
expression for w{r) within the hard-sphere core r < a, 
than the DFT used to obtain the results in Fig. [7] - see 
the discussion in Sec. |TT|). The boundaries of the two 
phase region lie either side of the A-line. As the value 
of e^^ is increased these coexistence lines meet at the 
A-line. This meeting point is a tricritical point. Above 
the tricritical points, the phase transition is second order 
and follows the A-line exactly. 

In Fig. [8] we display plots of the amplitude as a func- 
tion of p for a number of different values of e^^. For 
< e^^, we see that the amplitude varies discontin- 
uously at the phase transition, jumping from ^ = in 
the uniform phases to a non-zero value in the modulated 
phase. However, for > e^^ we see that the ampli- 
tude changes continuously across the phase transition. 
In Fig. [HI we plot the amplitude squared as a function of 
p for case when e^^ = 0.76 > e^^. We see that when 
A is small, (x {p ~ p\), indicating that the result in 
Eq. ([38]) , obtained from the order parameter theory, also 
applies to solutions of the full DFT. 
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FIG. 8: Variation of the order parameter, the amplitude A, 
as a function of the chemical potential /i, for various values 
of . The curves for = 0.75 and 0.76 are both con- 
tinuous and correspond to values of above the tricritical 
points. The curves for = 0.66 and 0.69 have discontinu- 
ities (marked by the dashed lines) at the phase transitions and 
these curves correspond to values of below the tricritical 
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FIG. 9: Plot of the amplitude squared ^ versus density p for 
the parameters A = 0.5 and — 0.76. The straight line is 
the line A? = 1.7(p - px), where p^a'-^ = 0.1835. 



VI. DISCUSSION AND CONCLUSIONS 

In this paper we have used DFT to study systems 
composed of particles with competing interactions. This 
model fluid exhibits phase transitions from the uniform 
fluid to a modulated fluid phase. In order to elucidate 
the precise nature of the phase behaviour of this sys- 
tem we developed an order parameter (Landau) theory 
for the phase transition, using the amplitude A of the 
modulations as the order parameter. Both the full DFT 
theory and the order parameter theory predict the fol- 
lowing phase behaviour: when the parameter e"""^ < e^^. 



on increasing the fluid density one finds that there is a 
first order phase transition from the uniform fluid phase 
to the modulated phase. At the transition point, there 
is discontinuous change in both density p and the am- 
plitude A, which jumps from ^ = to ,4 ^ discon- 
tinuously. However, when e^^ > e^^ , on increasing p 
one flnds that there is a second order (continuous) phase 
transition from the uniform fluid phase to the modulated 
phase. At the transition point, which is the A-line, both 
density p and the amplitude A vary continuously. Both 
outside and on the A-line itself, A ~ 0. On moving off 
the A-line, one finds that the amplitude A increases in 
a continuous manner. The precise form of this increase, 
for small values of A, is given by Eq. 

The phase behaviour predicted by the theory in this 
paper is in qualitative agreement with the results ob- 
served recently in Monte Carlo computer simulations of 
the present model p!3 |. However, one should bear in mind 
that the present theory is a mean field theory. The phase 
behaviour of the present system is somewhat analo gou s 
to that observed in diblock copolymer systems [l|, H l45j| . 
In these systems, the result is that when fluctuations be- 
yond mean fleld are taken into account, the second order 
transition becomes weakly first order [45|. This scenario 
is not supported for the present system by the Monte 
Carlo simulations in Ref. [l2|, but it may be the case 
that for some choices of parameters not explored in Ref. 

, that the transition at the A-line becomes weakly first 
order. 

We should also remind the reader that in the present 
study we have assumed throughout that the density pro- 
file varied only in one spatial dimension. Whilst this is 
true for any lamellar phase, one should expect the present 
system to exhibit cluster, bubble and perhaps other mod- 
ulated phases [l2|- The density profile for these phase 
will vary in more than one Cartesian direction. We plan 
to study the phases exhibited by the present DFT when 
the theory is not constrained to exhibiting modulations 
in only one direction. We expect the region in phase 
diagrams [H] and [7] labelled 'modulated fiuid' to be fur- 
ther subdivided into a number of different modulated 
phases, in a manner somewhat analogous to the mod- 
ulated phases displayed by a two dimensional fluid with 
competing interactions [2l|, [22] ■ 

The technological applications of fluids exhibiting 
modulated phases could be significant. For example, 
in display technologies or in making masks for micro- 
lithography. It is therefore important to understand and 
control the formation of the various type of modulated 
structures exhibited by these systems. The present study 
goes some way towards this goal. 
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APPENDIX: COMMON TANGENT 
CONSTRUCTION 

In this Appendix we show that the common tangent 
construction on the Helmholtz free energy per particle 
a = T/N yields the coexisting state points. 

From the free energy we may obtain the following 
quantities: (i) the pressure 



P 



dV 



N,T 



da 
dv 



(A.l) 



where v — V/N^l/pis the volume per particle and (ii) 



the chemical potential 



V,T 



da 
dv 



(A.2) 



If vi and V2 are the coexisting volumes per particle, then 
the conditions for mechanical equilibrium P{vi) — P{v2) 
and of chemical equilibrium n(vi) = Ijl(v2) give us 



da 
dv 



da 
dv 



a{vi) - a{v2) 

Vl — V2 



(A.3) 



which geometrically corresponds to the common tangent 
construction on a{v). 
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